
#include <iostream>
#include <sstream>
#include <vector>
#include <cassert>
#include <sys/time.h>

#include "geomTypes.h"

typedef float                                   Number;
typedef Point<3, Number>                        Point3;
typedef Vector<3, Number>                       Vector3;

#include "mlsComparison.h"

void createMesh(vector<Point3>&, vector<Point3>&, uint, uint);
void createControlPoints(vector<Point3>&, vector<Point3>&, uint);
double get_clock();
void printInfoAndExit();

int main(int argc, char **argv) {
   
   vector<Point3> P, Q;
   vector<Point3> positions, positionsIni;
   
   if (argc !=5) printInfoAndExit();
   
   unsigned int meshWidth, meshHeight, NCTRLPOINTS, method;
   meshWidth = atoi(argv[1]);
   meshHeight = atoi(argv[2]);
   NCTRLPOINTS = atoi(argv[3]);
   method = atoi(argv[4]);
   
   positions.resize(meshWidth*meshHeight);
   positionsIni.resize(meshWidth*meshHeight);
   
   createMesh(positions, positionsIni, meshWidth, meshHeight);
   createControlPoints(P, Q, NCTRLPOINTS);
   
   double t_start = get_clock();
   float p[3], pd[3];
   if (method==0) {
      for (unsigned int i=0; i<positionsIni.size(); ++i) {
         p[0] = positionsIni[i][0];
         p[1] = positionsIni[i][1];
         p[2] = positionsIni[i][2];
         mlsEval(P, Q, p, pd);
         positions[i] = Point3(pd[0], pd[1], pd[2]);
      }
   }
   else if (method==1) {
      for (unsigned int i=0; i<positionsIni.size(); ++i) {
         p[0] = positionsIni[i][0];
         p[1] = positionsIni[i][1];
         p[2] = positionsIni[i][2];
         mlsEval_OM_Horn(P, Q, p, pd);
         positions[i] = Point3(pd[0], pd[1], pd[2]);
      }
   }
   else if (method==2) {
      for (unsigned int i=0; i<positionsIni.size(); ++i) {
         p[0] = positionsIni[i][0];
         p[1] = positionsIni[i][1];
         p[2] = positionsIni[i][2];
         mlsEval_Quat(P, Q, p, pd);
         positions[i] = Point3(pd[0], pd[1], pd[2]);
      }
   }
   else if (method==3) {
      for (unsigned int i=0; i<positionsIni.size(); ++i) {
         p[0] = positionsIni[i][0];
         p[1] = positionsIni[i][1];
         p[2] = positionsIni[i][2];
         mlsEval_SVD(P, Q, p, pd);
         positions[i] = Point3(pd[0], pd[1], pd[2]);
      }
   }
   else if (method==4) {
      for (unsigned int i=0; i<positionsIni.size(); ++i) {
         p[0] = positionsIni[i][0];
         p[1] = positionsIni[i][1];
         p[2] = positionsIni[i][2];
         mlsEval_nearly(P, Q, p, pd);
         positions[i] = Point3(pd[0], pd[1], pd[2]);
      }
   }
   else cout<<"attention!!! method: invalid option!: "<<method<<endl;
   
   double mlsTime = get_clock() - t_start;
   //~ cout<<"  mlsTime: "<<mlsTime<<" sec."<<endl;
   cout<<mlsTime<<endl;
   
   return 0;
}

void createMesh(vector<Point3> &pos, vector<Point3> &posIni, uint mesh_height, uint mesh_width) {
   
   assert(pos.size() == posIni.size());
   for (unsigned i=0; i<pos.size(); ++i) {
      const float x = rand()/(float)RAND_MAX;
      const float y = rand()/(float)RAND_MAX;
      const float z = rand()/(float)RAND_MAX;
      pos[i] = Point3(x,y,z);
   }
   posIni = pos;
}

void createControlPoints(vector<Point3> &P, vector<Point3> &Q, uint NCTRLPOINTS) {
   
   srand(time(NULL));
   for (unsigned i=0; i<NCTRLPOINTS; ++i) {
      {
      const float x = rand()/(float)RAND_MAX;
      const float y = rand()/(float)RAND_MAX;
      const float z = rand()/(float)RAND_MAX;
      P.push_back(Point3(x,y,z));
      }{
      const float x = rand()/(float)RAND_MAX;
      const float y = rand()/(float)RAND_MAX;
      const float z = rand()/(float)RAND_MAX;
      Q.push_back(Point3(x,y,z));
      }
   }
}

double get_clock() {
   
   struct timeval t;
   struct timezone tz;
   gettimeofday(&t,&tz);
   return (double) t.tv_sec + (double) t.tv_usec*1E-6;
}

void printInfoAndExit() {
   
   cout<<endl;
   cout<<"Usage: mlsComparison meshWidth[0..n] meshHeight[0..n] ncontrolPoints[0..n] method[0,1,2] {nthreads[0..n]}"<<endl;
   exit(0);
}
